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Abstract 

We present local and global techniques to visualize three-dimensional vector field data. Us- 
ing the Line Integral Convolution (LIC) method to image the global vector field, our new 
algorithm allows the user to introduce colored “dye” into the vector field to highlight local 
flow features. A fast algorithm is proposed that quickly recomputes the dyed LIC images. 
In addition, we introduce volume rendering methods that can map the LIC texture on any 
contour surface and/or translucent region defined by additional scalar quantities, and can 
follow the advection of colored dye throughout the volume. 


+ This research was supported in part by the National Aeronautics and Space Administration under NASA 
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1 Introduction 


Visualizing vector field data is challenging because it is intrinsically difficult to visually con- 
vey large amounts of three-dimensional directional information. In fluid flow experiments, 
external materials such as dye, hydrogen bubbles, or heat energy are injected into the flow. 
The advection of these external materials can create stream lines, streak lines, or path lines 
to highlight the flow patterns. Analogies to these experimental techniques have been adopted 
by scientific visualization researchers. Numerical methods and three-dimensional computer 
graphics techniques have been used to produce graphical icons such as arrows, motion par- 
ticles, stream lines, stream ribbons, and stream tubes that act as three-dimensional depth 
cues. While these techniques are effective in revealing the flow field’s local features, the 
inherent two-dimensional display of the computer screen and its limited spatial resolution 
restrict the n um ber of graphical icons that can be displayed at one time. In addition, inter- 
active guidance that shows the user where to probe locally to obtain specific information is 
usually not available. 

Additional techniques for vector field visualization include global imaging techniques. 
Crawfis and Max [1] [2] proposed direct volume rendering methods to create images of entire 
vector fields. Vector kernels and texture splats are used to construct three-dimensional scalar 
signals from the vector data, van Wijk [3] proposed a Spot Noise method using stretched 
ellipses to create two-dimensional textures that can be mapped onto parametric surfaces. 
Max et al.[A) further utilized the spot noise method to visualize three-dimensional velocity 
fields near contour surfaces. Cabral and Leedom [5] presented a Line Integral Convolution 
(LIC) method, which uses a one-dimensional low pass filter to convolve a white noise texture 
based on the directional information of the vector field. These methods can successfully 
illustrate the global behavior of vector fields; however, little or no user probing capability is 
provided, so specific information about the local behavior of the field is limited. Furthermore, 
some of these methods are difficult to apply to unstructured meshes. 

In this paper, we present methods that integrate local and global visualization techniques 
to explore three-dimensional vector field data on regular grids. Using the Line Integral 
Convolution method as the underlying algorithm, we enable local probing by allowing the 
user to introduce “dyes” of various colors into the 2D/3D LIC flow field. The inserted dye 
propagates through the flow field, highlighting local flow features such as wavefronts, while 
the standard LIC texture still illustrates the global motion. We use three-dimensional direct 
volume rendering techniques that can map the LIC textures onto any contour surfaces or 
translucent region, and also can display the propagation of the dye through the volume. 

We begin the paper by giving an overview of the LIC algorithm. Next, we describe the 
method of dye insertion. We propose a fast algorithm that can significantly reduce the time 
needed to recompute the LIC images upon insertion of new dye material. We then describe 
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the convolution kernel and how to control the advection distance of the dye. We also describe 
direct volume rendering methods that are used to render three-dimensional LIC data. We 
conclude the paper by presenting results of applying the techniques on a variety of large-scale 
scientific data sets. 

2 Background 

In this section, we give an overview of the LIC algorithm originally proposed by Cabral and 
Leedom [5] . Although the method can be applied to both two- and three-dimensional vector 
data, in this section we restrict our discussion to the two-dimensional case to illustrate the 
basic algorithm. We omit the description of some recent improvements to the LIC algorithm 
subsequently proposed by Forssell and Cohen [6] and by Stalling and Hege [7]. However, we 
note that these newer methods can be easily adapted to our new algorithm. 

The LIC algorithm takes a vector field and a texture as inputs. The input texture is 
usually white noise data with the same resolution as the vector field. The output of the 
algorithm is a scalar field resulting from a. local blurring of the input noise texture. The 
LIC algorithm carries out the local blurring by applying a one-dimensional low pass filter 
convolution throughout the input texture. The convolution kernel follows the direction of 
the streamlines originating from the corresponding grid points of the vector field in both 
positive and negative directions. As a result, the intensity values of the output scalar field 
are strongly related to the vector field’s local flow direction. This convolution can be ex- 
pressed as follow: 

„ , , ELoCntPOAi + ElLoCnfCOk; 

FU " y) ~ eLoP+eLo'*; 

rsi+Asi 

hi = / K,{w)dw 

J Si 

Where 

• F out {x,y) is the output pixel value at point ( x,y ). 

• Fi n (Pi) is the input pixel value at point P,. 

• / and /' are the convolution distances along the positive and negative directions, re- 
spectively. 

• Pi represents the ith cell the streamline steps in the positive direction, and P( represents 
a step in the negative direction. 

• Po = (*,y). 


(2) 
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• hi and h' t are the weighting variables computed from the exact integral of the convolu- 
tion kernel k(w). 


• k(w) is the low pass filter used for the LIC. 

• A S{ is the arc length between the point s; and s; +1 along the streamline. 

• s 0 = 0. 

The low pass filter used in [5] is a Hanning ripple function, which has a period of 2 tt. By 
shifting the phase of the filter function while performing the convolution, the algorithm can 
generate a sequence of LIC images to create a periodic motion effect. 

3 Dye Injection 

The intensity value of each LIC output cell can be represented as an average of input texture 
values along the streamline. This leads to the result that each output cell has an intensity 
correlated with cells along the streamline, but not with other cells. The use of the white 
noise input to the LIC algorithm assures that the boundaries of neighboring streamlines are 
not obscured. In addition, the LIC algorithm changes the phase of the convolution filter to 
“push” the noise texture along streamlines. 

Through enhancing the individual streamlines and moving the noise texture, the LIC 
algorithm provides an excellent visual representation of the flow motion in the vector field. 
However, because of the nature of the LIC method, the correspondence between neighboring 
streamlines is difficult to observe. Sometimes, the user needs an observable correspondence in 
order to track local flow features, such as wavefronts. To resolve this discrepancy, we propose 
to use local flow visualization techniques, namely the introduction of foreign materials into 
the flow, and to observe the advection of these materials within the flow field. In the next 
section, we present a technique for injecting dye into the LIC field to highlight the flow field’s 
local features. 

3.1 Dye Smearing 

The LIC algorithm can be used to create a “smearing” effect by convolving the filter with an 
input image. This allows the user to generate motion blur and additional artistic effects, as 
demonstrated in [5]. The dye insertion method utilizes LIC’s natural “smearing” to simulate 
advection of dye within the flow field. We simulate the dye injection by assigning colors to 
isolated local regions in the input white noise texture. Cells whose streamlines pass through 
such regions receive color contributions from the dye. In addition, the phase shifting of the 
LIC algorithm can push the concentration of the dye along the streamlines. This creates the 
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effect of dye propagation. In the standard LIC algorithm, the convolution for each cell is 
performed in both positive and negative streamline directions to conserve symmetry. If we 
directly apply the standard LIC method to the input dyed texture, the upstream of the dyed 
area will be colored because the cells’ positive streamlines in that area will pass through the 
dyed area and obtain the color contribution. However, to create a correct motion effect, the 
dye should smear only in the forward direction of the flow field. Therefore the dye should 
color only those cells in the downstream direction that correspond to cells whose negative 
streamlines pass through the dyed areas. To overcome this difficulty, initially a regular LIC 
image F out (x , y) for each animation step is computed using the white noise input. When the 
user injects the dye, we apply the LIC convolution using the dyed texture input to those 
cells that will be affected by the dye. The convolution is applied along these cells’ backward 
streamline directions to ensure that the dye will only smear forwards. The results are then 
stored into D out (x,y). The final image of the LIC with dye advection can be obtained by 
using the formula: 


Final out (x , y) - D out {x , y) ® F out (x , y) (3) 

The operator 0 overwrites the standard LIC pixel values by the dyed values. 

An important issue that must be addressed stems from the fact that every time the user 
injects a new dye, only a small portion of the cells are affected by the dye. In the next 
section, we propose a fast searching algorithm to rapidly locate those cells affected by the 
injection of new dye. 

3.2 Fast Searching Algorithm 

From our previous description, we know that only cells whose negative streamlines pass 
through the dyed region within a distance of the convolution length receive color contribu- 
tions from the dye. We define a relation: 

Definition 1 For cells a and (3, a-^/3 if and only if there is a backward streamline flowing 
through a and stepping into f3 

Given a dyed cell a , the set of cells, 0(a), that are affected by this dye can then be expressed 
as: 

Axiom 1 f 1(a) = {p\p-¥a : d(p, a) < L} 

where d(p, a) is the arc length of the streamline from point p to point a, and L is the 
convolution length. 

Given a cell a, the goal is to quickly locate the set of cells 0(a) without searching through 
the entire field. A brute force method is to create a list for each cell. The list contains all the 
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cells that step backwards along streamlines into that cell. This information can be gathered 
while we first perform the standard LIC, since all the necessary streamlines are computed at 
that time. To decide upon the set ft(a), we can simply traverse through the list and retrieve 
the cells. 

While the above method is simple, the memory requirement to perform it would be 
overwhelming. This is because each cell's list would contain many cells along the streamlines. 
To resolve this, we have designed a new algorithm. First, we define a new relation: 

Definition 2 For cells a and /3. a A (3 if a—>(3 or 3{po,pi, ■ ■ ■ p n }, such that 

a-±po,po-*pu • • • pn~>(3- 

From definition 1 and definition 2, we know that: 

Axiom 2 if a-¥0 then a A /?. 

Note that the above axiom does not necessarily hold in the reverse direction. Secondly, 
define a set u) as: 

Definition 3 u(a) = {p\p A a, d(p, a) < L] 

From axiom 2, we know that: 

Axiom 3 w(a) D 0(a). 

Our fast searching algorithm can locate the cells in the set u(a) that are a superset of 0(a). 
In practice, the difference between the number of cells in the set 0(a) and the set w(a) is 
very small. To locate the cells in cu, we define a new relation: 

Definition 4 A cell a -4 (3 if and only if a-»/? and a is /3 ’s direct neighbor. 

We call this relation a direct flow-back relation. The direct neighbors of a cell are those 
cells that are directly adjacent to that cell in the physical space. For example, a cell has 
eight direct neighbors in a two-dimensional Cartesian grid, and has 26 direct neighbors in 
three-dimensional space. 

The direct flow-back neighbors for each cell can be obtained and stored when we first use 
the standard LIC algorithm. Given a root cell a, the cell and its recursive direct flow-back 
neighbors actually constitute a directed graph. The nodes in this directed graph are then 
the members of the set u>. Therefore, the u>(a) can be found by using: 

Algorithm 1 A member of tu(a) can be found by using a Breadth First Search method 
starting from cell a and traversing through its direct flow-back neighbors recursively. 
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Figure 1: Streamlines and dye 

Figure 1 shows an example of a. flow field. The shaded cell is the original dyed cell. 
Figure 2 illustrates the resultant directed graph if we only consider the streamline A, B , and 

C. 

Our new algorithm is considerably more memory efficient than the brute force method 
because only the direct flow-back neighbors need to be stored for each cell. The number of 
direct flow-back neighbors will usually be much lower than the actual number of the direct 
neighbors for a. cell, with the exception of cells near critical points. We have found that 
the average number of direct flow-back neighbors for each cell in our experimental three- 
dimensional flows is approximately two or three. 

3.3 Convolution Kernel 

The Hanning windowed filter used in [5] can not be easily adopted by our algorithm to create 
natural dye advection effect. This is because that the Hanning ripple function is a periodic 
function with multiple peaks, and the Hanning window function attenuates the values of 
the incoming function values at two ends. Even though several parameters of these two 
functions could be modulated, it is still difficult to construct a good combination of the 
parameters such that the concentration of the dye can advect along the streamline smoothly. 
In our implementation, we use a simple box filter instead of the Hanning filter so that the 
appearance of the dye advection can be controlled more easily. The equation 2 becomes: 

rSi+Asi 

hi — / 1 dw = A Si (4) 

To add the ability to shift the filter box to create flow motion in the animation sequence, 
we revise the equation 1 to: 
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Figure 2: Flow-back directed graph 


F 0 ut(l,x,y) 


ELo c,4C)/».r(/, »■) + E'lo jUjM r ( j , 

EU'>.r(/,si) + ELoASr(/,si) 


Where 


T(I,Si) 


0 otherwise 


• $ is the number of LIC animation steps 

• /€[0,$-l] is the current animation step 

• A L is the length of the filter box 


(5) 


Figure 3 shows the phase shift of the box filter. 

In equation 5, T(/, s*) is a predicate that determines if the current pixel is within the 
range of the filter box L , £ L + A L). Note that if + A L is greater than L, this value 
needs to be wrapped around so that it become L + AL) mod L. The new notation G* 
that we use above is therefore based on this relationship. 

It is known that the box filter could cause artifacts when the boxes reenter the interval. 
A solution to avoid this problem can be found in [7]. 


3.4 Advection Length 

In the standard LIC algorithm, there is only a single convolution length defined for all cells 
in the field. 1 This global convolution length determines the distance that the dye can travel 

*To simplify our explanation here, we omit the variation that the convolution length can be scaled based 
on the vector magnitudes. 
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and this distance is often desired to be long enough so that the local flow features can be 
clearly highlighted. However, extending the convolution length globally would decrease the 
contrast of the LIC texture and slow down the computation. In our algorithm, we allow the 
user to specify convolution lengths separately for the dye advection. 

To use a longer convolution length for the dye advection, more LIC animation steps are 
needed to complete the animation. Assuming that the global convolution length is L, the 
number of animation steps for the standard LIC convolution is $ and that the convolution 
length for the dye is T, then X = $ x £ is the number of animation steps required. Here we 
assume that T is always a multiple of L. We revise the LIC formula for the dyed advection: 
D out (I,x,y ) = 

Z-loDa^hiT'iLsi) 

ErJo hiT’{I,Si) 


Where 




i if Si e* [|r,|r + AL] 

0 otherwise 


• Di n (Pi) is the input dyed texture 

• D 0Ut (I,x,y ) is the output value 
for the dyed cell(;r, y) at step I. 

• / €■ [0,Z — 1] is the current animation step 

• A L is the length of the filter box. 


Note that for those cells that are not affected by the dye, only $ steps of LIC output 
are computed. To combine these standard LIC outputs with the results from the LIC dye 
advection, we can use the formula: 


8 



( 7 ) 


Final 0Ut (I , x, y) = D out {I , x , y) <g> F out (/ mod $, x, y ) 
where the meaning of the operator <g> is explained in section 3.1. 

4 3D LIC Rendering 

While the output of the three-dimensional LIC algorithm is a scalar field, standard visu- 
alization techniques, such as isosurface extraction or volume rendering, cannot be directly 
applied to the output. This is because the values of the LIC volume consists primarily of 
noise signals that are used to construct the vector textures and do not themselves have any 
physical meaning. In this section, we present direct volume rendering methods that can 
render three-dimensional LIC textures on arbitrary surfaces or regions derived by additional 
scalar quantities, and also can display the propagation of the dye through the volume. In 
the sections that follow, we first introduce a bi-variate volume rendering model and then 
describe the methods used to image the dye propagation. 

4.1 Bi- Variate Volume Rendering 

To render the LIC texture onto local regions, our volume rendering process takes two volume 
data sets with separate transfer functions as the input. A scalar variable, such as the 
magnitude of velocity, or pressure, or of vorticity, is used as the primary data. The opacity 
map of the primary data, which is manipulated in the same way as in the standard volume 
rendering method, is used to define the transparency of each voxel in the volume. The LIC 
texture, which serves as the secondary data , has its own color map but uses the transparency 
defined by the primary data. To render these two input data sets, we define a volume mixture , 
which is similar to what is used in [8], to blend these two input data sets. The color of each 
voxel V, Cvi comes from the contributions of both data sets and can be expressed as: 

Cv = 0xC p + {l-0) xd s , (8) 

where C p is the color from the primary data at voxel V and C s is the color from the secondary 
data. The (3 E [0, 1] is a weight between the primary and secondary data specified by the 
user. When /3 is 1, the volume mixture is equivalent to the scalar volume; when (3 is 0, the 
volume mixture is the same as the LIC volume. A gradual shifting of this parameter from 
1 to 0 enables the user to track features from the scalar data at the beginning and visualize 
the counterpart of the vector data and vice versa. Our new rendering model allows the user 
to interactively change the weighting parameter and visualize the results. 
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Our bi-variate volume rendering method uses the standard front-to-back compositing 
method, where the final color at pixel p can be expressed by the following formula 

C, = ±W(i)Cv„ (9) 

2—1 

where Cy, is the color contribution at the i-th sample point Vi. W(i) is the light attenuation 
factor, which is computed from the formula: 

W(i) = a(V z )[l-Y / W(j)}. (10) 

3 - 0 

o(V)) is the opacity at the point V t and is defined by the primary data set’s opacity map and 
its data value at that point. We can replace the CV- in equation 9 with the volume mixture 
in equation 8 and obtain: 

C„ = t. w (>)Cv, = 

i— 1 

•£ W(i)\p X Cv M + (1 - ft) X CV iW ] = 

2—1 

+ (1 - l3){£,W(i)Cv M }. (11) 

i= 1 t= 1 

Here Cvi( v ) is the color contribution at point V } from the primary data set and CV ; ( S ) is 
the color contribution from the secondary data set. Equation 11 reveals a very interesting 
property: when we render the volume mixture, we can render each volume separately. The 
mixture image can be generated by a simple two-dimensional image blending. 

Thus the bi- variate volume rendering model can be used to intermix the scalar and vector 
information. Moreover, the user can interactively adjust the volume blending parameter f3 
to track important features from one variable and visualize the counterparts of the other 
variable. 

4.2 Dye Rendering 

It is helpful to display the LIC texture onto local regions such as isosurfaces and translucent 
areas. However, this could give the user the wrong impression that the flow is always along 
the surface or local regions. This problem can be remedied if we can visualize the dye 
propagation through the volume. In our method, we allow the user to specify opacity of the 
injected dye. The opacity values in the LIC dye input are convolved in the same way that we 
convolved the pixel values. As a result, the output volumes of the LIC computation consist 
of both color and opacity information. 
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To take the dye opacity into consideration when performing volume rendering, the bi- 
variate rendering model needs to be slightly modified. This is because the LIC volume, as 
a secondary data, also affects the volume’s opacities. When rendering the LIC volume, we 
calculate the opacity value at each sample point as the sum of the opacities from the primary 
data and from the LIC volume. The color contribution of the LIC voxel at that sample point 
is then computed using this opacity. Following is pseudo code for rendering the LIC volume 
using ray casting. 

for (each ray from the image plane) { 
for (each sample step along the ray) { 
opacity = primary opacity + dye opacity; 
accumulate the LIC color; 
accumulate the opacity; 

} 

> 


A simple but useful extension is that the color of dye can be determined by additional 
scalar variables. In this way both the flow direction and scalar distribution can be visualized 
at the same time. 

5 Results and Discussion 

We have implemented our algorithms using C++ and OpenGL. The software combines 
modules of Line Integral Convolution with Dye Injection and Bi-variate volume rendering. 
In our current software implementation, the user can inject dye by slicing through the three- 
dimensional LIC output and specifying the dye locations. Our fast search algorithm can 
rapidly recompute the LIC output based on the dye position and then feed the results back 
to the volume Tenderer module. We have found that our fast search algorithm can reduce 
the compute time, on average, by over 90% when compared to the standard LIC method 
for recomputing the animation sequence of the LIC data. The output of the LIC volume 
contains red, green, blue and opacity components for each voxel. For the volume rendering 
module, we define a C++ volume Tenderer class. The Tenderer instance of this class can be 
either a primary data Tenderer or a secondary data Tenderer. Through a software pointer 
between the primary data Tenderer and secondary data Tenderer, the secondary Tenderer 
can access the primary Tenderer’s transfer function, volume data, and other auxiliary data 
structures to facilitate the bi-variate rendering model. A secondary Tenderer class instance 
is created for each LIC output in the animation sequence. These Tenderer objects can be 
dispatched to different processors in a multiprocessor environment to perform a coarse grain 
parallel volume rendering. 
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Figures 4-6 show several LIC images with dye propagation. Figure 4 shows three different 
colored dye propagating through a two dimensional vector field. Figure 5 is an animation 
sequence of dye advection in a 32 x 30 x 30 vector field from a three-dimensional combustion 
simulation. Figure 6 is an animation sequence of dye advection in a 96 x 96 x 96 tornado 
data set. The dye is colored based on the magnitude of the vector field. 

6 Conclusion and Future Work 

We have presented new algorithms to visualize three-dimensional vector fields using Line 
Integral Convolution methods. Our algorithm allows the user to insert dye into the flow field 
to enhance local flow features. Given a region with dye, our algorithm can rapidly locate 
the cells that are on the path of the dye propagation. Therefore, the LIC computation needs 
to be applied only on a small subset of the entire field. We proposed a bi-variate volume 
rendering method to display three-dimensional LIC textures on arbitrary contour surfaces 
or translucent regions that are derived from additional scalar quantities. The new rendering 
methods can also be used to render the dye propagation through the volume. 

Future work includes providing more flexible control for modeling the behavior of the dye. 
This includes using different filter kernels to control the appearance of the dye. An additional 
useful feature would be to provide a natural and effective way for the user to probe the three- 
dimensional space and inject the dye. Finally, we plan to utilize the coherence between LIC 
volumes in an animation sequence to speed-up the volume rendering and allow for interactive 
data exploration. 
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2D LIC images with dye 


Figure 6 : Dye advection in a 3D tornado simulation 


Figure 4 


Figure 5: Dye advection in a 3D combustion simulation 
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